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The Landau-Lifshitz equation for the magnetization dynamics of a single-domain magnetic system 
is solved using the methods of self-organization. The description takes into account the torque due to 
spin transfer. The potential energy of the system includes the uniaxial and easy-plane anisotropics, 
and the Zeeman energy due to an external magnetic field. The equilibrium and stationary states 
are investigated as a function of the spin current and external magnetic field. The presented bifur- 
cation diagram allows to determine the margins of a neutral stability mode of the equilibrium and 
stationary states for different values of the easy-plane anisotropy constant. Using the power spectral 
density method, the trajectory tracing, Hausdorff dimension, and maximum Lyapunov exponent, 
the dynamics of the phase states in an external magnetic field is demonstrated. The analytical 
transcendent equations for switching between different equilibrium states are also obtained, proving 
the importance of phase averaging. 
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I. INTRODUCTION 

Spintronics is a new and intensively developing field of 
physics, which attracts significant scientific interest from 
both fund amen tal an d application points of view (see, 
e.g., Refs. |ll2l l3ll4l l5l6| and the references therein). The 
main objective of spintronics is search for new materi- 
als for devices of future generation (like, for instance, 
spin transistors, quantum computing devices, solar cells, 
etc.) and with fundamentally new operation princi- 
ples in comparison to traditional devices based on the 
electron charge transport^"— One of the main prob- 
lems, which still require solution, concerns effective con- 
trol of spin polarization and spin injection in spintron- 
ics devices iSiiifi"— This problem triggered an exten- 
sive theoretical and experimental search for new ma- 
terials and devices, like metallic systems consisting of 
ferromagnetic and nonmagnetic layers or tunnel junc- 
tions including ferromagnetic and diluted magnetic (non- 
magnetic) semiconductors Significant progress has 
been made recently in the technology and physics of 
doped semiconductors with spontaneous spin polariza- 
tion, and in the spin injection from magnetic to non- 
magnetic systemsii^"— Several theoretical models have 
been also proposed to describe doping-induced ferromag- 
netism and spin injectioni^'^i^'*i^^i^^i^''i^^i^^i^° 

Theoretical studies carried out to date usually treated 
the spin subsystem by techniques developed for equi- 
librium situations. However, such an approach applied 
to the spin-transport phenomena ignores the fact that 



both spin polarization and spin injection are dynamical 
processes and therefore should be treated by methods 
of non-equilibrium thermodynamics^"— This is partic- 
ularly important for open dynamical systems that are 
able to exchange particles and energy with the environ- 
ment (which is true for the spin subsystem considered 
here). The dynamical processes may lead to stationary 
non-equilibrium states with high ordering level. Such co- 
operative phenomena in the system of magnetic ions and 
mobile band electrons can take place far from an equilib- 
rium state. Generally, the system can have several sta- 
tionary solutions described by the same set of nonlinear 
time-dependent differential equations, which are different 
regarding their stability against small fluctuations. 

In the case of stable solution, fluctuations are sup- 
pressed within a short period of time. If a solution is un- 
stable, fluctuations can be amplifled and the system may 
be switched to another stable state. The new state may 
have the same or even higher ordering level, when it is de- 
scribed with a lower symmetry group. Such spontaneous 
switching between different possible states is known as a 
self-organization of the system. The methodology based 
on non-equilibrium thermodynamics, developed for in- 
vestigating self-organization processes, offers promising 
approach to such spintronics problems like spin polariza- 
tion and spin injection. This may lead to new quantita- 
tive and even qualitative results, which are not accessible 
within the traditional approaches. 

The problem studied in this paper concerns the dy- 
namics of a single-domain magnetic system interacting 
with a spin current flowing through it - the problem al- 
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ready investigated in detail by Siina* using the standard 
methods. The approach based on the self-organization 
techniques yields some new interesting results. 

The paper is organized as follows. The model is de- 
scribed in Section 2, where the symmetry of a quasi-one- 
dimensional ferromagneli^ is analyzed in detail. Analyt- 
ical solutions corresponding to four equilibrium and two 
independent stationary states are also derived there. The 
stability analysis regarding the influence of small pertur- 
bations shows that the spatial angles describing the mag- 
netization vector can vary either monotonically or in an 
oscillatory manner with time, depending in a complex 
way on the system parameters. 

Section 3 presents the results of our numerical calcula- 
tions and the relevant discussion. This section is split 
into two subsections in order to separate the analysis 
of equilibrium and stationary states from the analysis 
of the time evolution of the system. We also discuss 
the magnetization switching and its dependence on the 
magnetic field and spin current. It is shown that the 
equilibrium and stationary states are characterized by 
a different dynamics of the oscillatory modes which can 
be excited by small perturbations of the phase variables. 
The bifurcation points of the system (the points where 
the system changes its behavior qualitatively) are also 
determined, and their dependence on the external pa- 
rameters is investigated. The second subsection presents 
the results obtained from the numerical simulations. A 
significant attention is paid to the analysis of the sys- 
tem evolution through several types of different mag- 
netization precession modes. It is shown that the sys- 
tem can be switched between different equilibrium states, 
ruz = ±1, in a controllable way, which is of particular 
interest for the spintronics applications. The detailed 
analysis of the system evolution has been performed us- 
ing the self-organization methodology; namely, the den- 
sity histogram of the longitudinal magnetization compo- 
nent, power spectral density, Hausdorff dimension, maxi- 
mum Lyapunov exponent, and the phase trajectory trac- 
ing curve. This allowed us to clarify peculiarities of the 
magnetization orientation for different values of the ap- 
plied magnetic field. The derived analytical expressions 
for the magnetic switching between the states rriz = ±1 
made it possible to analyze the dependence of the switch- 
ing time on the external magnetic field and the spin cur- 
rent. The last Section 4 includes the summary and main 
conclusions. 



II. THEORETICAL MODEL AND METHOD 

We consider a model system studied recently by Sun.^" 
Accordingly, we assume a single-domain quasi-one di- 
mensional ferromagnet of length l„i along the axis x and 
a square cross-section (a x a) in the plane y — z. The sys- 
tem is described by the easy-axis (along the axis z) and 
easy-plane (y — z) magnetic anisotropics. In addition, an 
external magnetic field H is applied in the plane y — z 




FIG. 1: Model geometry and main variables of the system. 



at an angle ij) with respect to the axis z. The external 
field H and the easy-plane anisotropy are described in di- 
mensionless relative units, h — H/Hk and hp = Kp/K, 
where Hk = 2K/M, K and Kp are respectively the easy- 
axis and easy-plane anisotropy constants, and M is the 
absolute value of magnetization. 

The system interacts with a spin current flowing 
along the wire. The incident spin current is described 
by the two parameters, 77 and n^. The parameter ry de- 
scribes the degree of spin polarization of the incoming 
charge current /, whereas ng is a unit vector along the 
corresponding spin polarization. We assume that n^. is 
in the y — z plane, and forms an angle (j) with the axis z. 
As in Ref. 20], the spin current is described in dimen- 
sionless units as = nshrjl / AelmO^ K . In the following 
we restrict considerations to homogeneous dynamics, and 
assume that the system absorbs the perpendicular to the 
magnetization component of the incident spin current. 
This produces an additional torque acting on the system 
magnetization. The geometry of the system, together 
with the reference frame used in the theoretical descrip- 
tion, are shown in Fig. 1. 

As in Ref. 20], we assume that the energy den- 
sity U includes the energy of uniaxial magnetic 
anisotropy, Uk = K sin^ 9, the energy of easy- 
plane anisotropy. Up = Kp {sin^ 9 cos^ (p — 1), and 
the Zeeman energy due to an external magnetic field, 
Uh — —K{h± sin 9 sin tp + cos 9) , where h± — h sin tp 
and hz — h cos ?/; are the components of the vector h (see 
Fig. 1). Thus, one can write 

C/o(0,^)^^^^Z(^) sin^fl 

— 2 (/iz cos6' -f- /i^ sin6'siniy9), (1) 
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where 



Z{(p) = (1 + hpcos^ 



(2) 



and a constant term in Eq. (1) has been omitted. 

The homogeneous dynamics of the magnetization M is 
determined by the Landau-Lifshitz equations^^ includ- 
ing also the torque due to the spin transfer^Si 
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with 

Ai = Z{lp) cos9 + h 
1 



sin6'(a^i + A2), 
=aA2-Ai, 

OT 



(3) 



sm ( 
1 



- {h± cos 9 sin (p + hs± cos (p) , (4) 



A2 = —hp sin 2ip H {h± cos ip — hs± cos 9 sin ip) . 

2 sin 9 

Here, h± = hsiiiip and hz = hcosip are the applied 
field components, hs± = /is sin and hsz = /is cos are 
the components of spin current (see Fig. 1), whereas 
T = t/[(l + a'^)/^Hk] denotes the dimensionless time, 
which depends on the gyromagnetic ratio 7 = g^b/h, 
the anisotropy field Hk, and the damping coefficient a 
{a ^ 1). When writing down the set of Eqs. Q, we 
assumed that the magnitude of the magnetization vector 
is constant, and g = 2. 

It is rather difficult to obtain analytical solution of 
Eqs. lO and in a general case. Therefore, in this 
paper we consider a particular case, when the external 
magnetic field h and the spin current hs are along the 
axis z of the reference system shown in Fig. 1 (similarly to 
the situation considered in Ref. |20] ) . Thus, we put h± — 
0, hz = h, hs± = and h^z = /is- Note that from now on, 
h and /is denote the z-components of the magnetic field 
and spin-current. Accordingly, both h and hg can take 
either positive or negative values. Consequently, Eqs. ^ 
and ^ can be rewritten as 



Uo{9,ip) = Z((^)sin2 6'- 2/1COS6I, 



(5) 



and 
89 
8^ 



8ip 
8^ 



a[Z{ip) cos9 + h] 



sm 2ip 4- /is 



■ sin 2ip + /ij 



[Z{lp)cos9 + h]. (6) 



Equations © present a generalization of Eq. (11) from 
Ref. jl^l to arbitrary values of the angle 9 (in Ref. |23| 
only the case of 16*1 ^ 1 was considered). 

The system described by Eqs. Q is invariant with re- 
spect to the following substitutions 

+ ^, (7) 



and 



p) ^ -K + ip, 

h ^ h, /is — >■ —/is, 

9^TT-9, 
p TT — p. 



(8) 




Magnetic field h [d.u.; 



FIG. 2: Possible equilibrium states of the system (here and 
in the following d.u. means dimensionless units). 



The property lO allows us to reduce the interval of the 
variable p to < p < tt, while the relations JHJ can be 
used to verify analytical results, in particular to check 
the solutions of the characteristic equations. 



A. Equilibrium solutions 

The equilibrium orientation of the magnetic moment 
M in the absence of spin current can be found from the 
conditions^ 



dUo 
89 



— 2sin6'e [Z{pe) COS^e + h] = 0, 



dUo 
dip 



-hp sin 



, sin 2(y9e 



(9) 



^From now on we label the equilibrium solutions with the 
index " e" , and the stationary ones with the index " 0" . 

Assuming that hp > Or^we find the equilibrium solu- 
tions of Eqs. (0) 



(1) 
(2) 

(3) 



sin^e = 0, Pe — arbitrary, 

C0s9e = —h, (/3e = 7r/2, (10) 

cos = - ^ ^ , , Pe^O. 
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As follows from Eqs. (|10|) . in the absence of magnetic 
field, h — 0, the system has the following energy states 
(phases): spin-degenerate states m — (0,0, ±1) corre- 
sponding to the energy Uq = 0, the state m = (1,0,0) 
with the energy Uq = 1, and the state m = (0,1,0) 
with the energy Uq = 1 + hp. Here, m is the unit vec- 
tor along the magnetization, m = M/M. The exter- 
nal magnetic field removes the spin degeneracy, and the 
level Uq — splits into two sub-levels with the corre- 
sponding energies Uai = — 2/i for the state m = (0,0, 1) 
and U02 = 2/1 for the state m = (0,0,-1). The mag- 
netic field also leads to precession of the ruz compo- 
nent around the direction of magnetic field for the states 
m = (Vl — /i^, 0, —h) with energy U03 = 1 + h'^ and 

m = (0, -^1 - (/i/(l -I- hp)f,-h/{l + hp)) with energy 
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Uq4 = 1 + hp + + hp). The state of energy U03 

exists in a certain range of magnetic fields, < \h\ < 1, 
while the state Uq4 corresponds to < \h\ < 1 + hp. It 
is worth noting that the first solution in Eqs. (|f 0|) cor- 
responds to the minimum of interaction energy for any 
value of the magnetic field h, while the third solution cor- 
responds to the maximum of energy for < \h\ < 1 + hp 
(Fig. 2). 

According to Ref. [33, a particular phase can be con- 
sidered stable in a certain range of the field h, provided 
the following conditions are obeyed: 



sm 



dHJo 
80^ dip^ 



>0, 



80 dip 



> 0. 



(11) 



When one of these expressions turns to zero, with a con- 
sequent change of sign, the corresponding phase becomes 
unstable. Substitution of the solutions H1U|) into Eqs. 
allows one to determine the conditions of the system sta- 
bility loss (for certain values of (p) for the states, which 
are odd with respect to the magnetic field, i.e., for the 
states Uoi and Uo2- The phases U03 and U04 (even at 
nonzero h) are unstable in any case. This indicates on 
the possibility of spin re-orientation leading to transitions 
between the existing phases under applied magnetic fields 
< /i < 1 or f < < 1 -h/ip. 

The boundary corresponding to the transitions be- 
tween two phases is defined by the following equations:^^ 



BU 

= 2sin0[Z{Lp)cos0 + h] = 0, 

da 



2[Z(</7)cos26l-|-/icos6'] = 0. 



(12) 



As the single solution of Ea. (|12() coincides with the 
first solution in (|10|l . it means that = Q, ip = 
0.5 arccos[-2(/i + l)/ft.p-l] for - (l + /ip) < /i < -1 with 
Uqi — —2h and = tt, ip = 0.5 arccos [2 {h — l)/hp — 1] 
for 1 < ft, < 1 -I- /ip with Uq2 = 2h describe the lines of 
the second order phase transition;^ in accordance with 
the Landau theoryi^ 



B. Constant energy solutions 

Using Eqs. Q and we find that the variation of 
Uo{0, ip) with time obeys the following equation: 



dUo{0.p) 



^ -2 sin^ 0{a [{Z{p) cos + h)"' 



dr 

-|-(0.5/ipSin2.^ 4- /i^)^] -I- hs [Z{p) cos 6* 

+h- a{Q.bhpSin2p + hs)\} . (13) 

^From this follows that the potential energy of the system 
is an integral of motion [Uo{0,ip) = Uq = const] when 



both a = and hg = 0, which allows to transform Eq. © 
to the form 



dp 

a7 " ^ 



Z{p) COS 



UoZ{p) + Zip) 
dip 



21 1/2 



(14) 



The " and " — " signs on the right-hand-side correspond 
to the two possible solutions for cos0 from Eq. (O at 
Uq =const. 

Equations (|14|l have real solutions only when 
h^ + Z{p)^ > UoZ{p), 

1 , . . ^^j.^ 



T < cos ( 



Zip) 



1? + ^ 



< 1. 



The analysis of condition (|15|l shows that the solutions 
of Eqs. ini) exist for \2h\ < Uq < C/03 and C/03 < C^o < 
Uoi, when it can be expressed by the Jacoby elliptic 
functions."*^ However, if the interaction energy Uq coin- 
cides with the energy of one of the equilibrium states H10|l , 
the solutions can be simplified to the inverse trigonomet- 
ric or exponential functions. In particular, for the case 
of Uq = ?7oi,2 = ±2|ft|, the solutions of Eqs. are 



(/? = =p arctan 



cos — ± 



1 - 2\h\ 



C(l 
e^(l^ 



e/C tanh(VCeT) 

C + gtanh^(^/CeT) 
hp) 



hp — \h 



etanh^(VCer) 

C^\h\-1. 



(16) 



Here, the upper sign corresponds to the interval of 1 < 
h < I -\- hp, while the lower one is for the interval of 
— (1 + ftp) < h < —1. When deriving Eq. ((TE|l we assumed 
the integration constant tq = 0. Note that for \h\ 1, 
the values of p and cos tend to the finite limits 



cos6'|,i|^i = ± 



=p arctan [hpT), 
_1 + 2''" 



1 



hy^ 



(17) 



It is important to note, that apart from the states 
described by Eq. (|16|1 . there also exists the state with 
cos 6^ = —1 for arbitrary ip and ft. > 0, and the state 
with cos 6* = 1 for ft. < 0. As follows from H17(l . if we 
neglect the plane anisotropy (ftp = 0) or take the limit 
of r — > 00, Eq. Hlt)|) for cos 9 gives the independent on p 
expressions mentioned above. In other words, Eqs. I|16|) 
and 1)17(1 prove that along the lines corresponding to the 
phase transitions of the second type, there are two pos- 
sible solutions for the longitudinal magnetization com- 
ponent: a constant one, and a periodic solution with a 
complex time dependence. In particular, for |ft| > 1 it is 
characterized by the cyclic frequency 



C.(ftp,ft) = [(|ft|-l)(l-Hftp-|ft|)] 



1/2 



(18) 



It is worth noting that Eq. ((TB|l for (^(t) and Eq. l(TS|) 
almost coincide with similar expressions for </?(t) and fre- 
quency Ljjp from Ref. po| , if one transforms the field in- 
terval from ft > — 1 to 1 < |ft| < 1 -I- ftp considered here. 
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When Uq = /7o3,4 the solution of H14() has the form: 

a) for the case Uq = Uq^: 

sin^ (p = cot^ {k'hp{T - To) /k) - fc2] , (19) 
Z{ip) cos 01,2 = -h± ^{Zi^)-l){Z{^)-h^), 

where k"^ = hp/{l + hp ~ h"^) and fc'^ = 1 - fc^. 

b) for the case C/q = ^-^04 and h > ^/T+Ti^: 

^ M^tan^ [/ip(T - To) /A] + 1 

Z((^)cos6ii,2 = -/i± (1 + /ip)"^/^ (20) 
^^{1 + hp- Z{ip)){h^ - (1 + V)^(^)), 

where = /ip(l + Zip)/ [(1 + hp)'^ - h'^] > 1, and the 
expression for Z{ip) is given by Eq. The integration 
constant tq in Eqs. (|19() and (|20|l can be defined from 
the initial conditions. The comparison of formulas (|17|l . 
(|19|1 and H2U|I proves that the magnetization dynamics 
corresponding to the equilibrium energy is expressed with 
the complex functions of h and hp, depending on the form 
of solutions l(Tn|l . 



C. Analytical solutions in the presence of spin 
current 

Let us now turn back to the discussion of the effect of 
spin current. It is quite natural to expect that the spin 
current can affect the system in a way similar to that of 
an external magnetic field. In a general case, the influ- 
ence of spin current on the magnetization dynamics can 
be studied only by the numerical integration of Eqs. JH)). 
This will be presented later. Below we consider three 
particular cases oi hg 0, for which some analytical so- 
lutions can be derived. 

The simplest case is the one with a = and h = 0, 
when Eqs. I© transform to 

— = - f - sin 2<^ -I- ft-s j sin 9, 



dip 
97 



-Z{tp) cost 



(21) 



Excluding the time variable, the above system of two 
equations can be reduced to a single equation with the 
corresponding solution 



= arcsin [Uo{ip)/Ziip)]^^\ 



(22) 



where 



/ 2hs arctan(^l + hp cot (p) \ 
C/o(^)=exp(^ ^jYTTp )■ ^''^ 

The integration constant in Eq. (|22|l was determined from 
the condition that 9 = 7r/2 corresponds to (p = 7r/2. 



Substituting Eq. into Eqs. and ((SJ, we obtain 



dip 

a7 



-Z{p) 



Zip,) 



1/2 



(24) 



It is obvious that the solution of Eq. H24|l can be found 
only by numerical methods. However, some analytical 
results are achievable regarding the magnitude and sign 
of hs- Equation H24(l can have real solutions only when 



^/T+hp InZjip) 
2 arctan {^1 + hp cot ip) 



> 



(25) 



If hs > 0, the condition 1251) is automatically fulfilled. 
With increasing spin current, the variables ^(t), 9{t) and 
Uo{t) relax to their equilibrium values, corresponding to 
the limiting case of /ig — > c» : 



ip = — arctan 



Vl + hp cot (^y/l + hpT - 

sin 6* = 0, C/o = 



(26) 



0. 



The integration constant for (/3(r) in Eq. H2t)|) was chosen 
to keep (p{t = 0) = 0. It is important to note that 
the nature of the solution H26|l coincides with the first 
solution in Eq. ((Tn|l . 

For hs < and for certain values of spin current, 
Eq. (|25|l is no longer valid and the system cannot be 
integrated. When the left and right sides of Eq. H25|l are 
equal, the system exists in the state 



sin0 = 1, (p ^ 0, Uo = l, 



(27) 



which coincides with the second solution in Eq. (|10|) . 

Therefore, when no magnetic field is applied, the spin 
current can excite the magnetic system, which in turn 
may relax either to the first or to the second state de- 
scribed by Eqs. (|10|) . depending on the sign of hs- In 
other words, the behavior of the system depends on both 
direction and magnitude of the spin current. One has to 
pay special attention to the exponential character of this 
dependence, because a comparatively small change of hs 
in certain direction can result in a significant change of 
the magnetization orientation and of the interaction en- 
ergy. It is also worth noting that the dependence of the 
magnetization components and the interaction energy on 
the applied magnetic field is much weaker, being linear 
or quadratic in h. 

The second case, where some analytical solutions are 
possible, is the situation with hp — 0. Equations © can 
be then rewritten for = cos 9 as 



dr 



a{l~ml){m, + h + hs/a) (28) 
dp 



dr 



'{rriz + h — ahs) 



It can be shown that for the case ah -\- hs ^ ±a, the 
solution of Eq. 128|l has the following form 



1 + 



(to^ -I- py 
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1 + TO 



zO 



,-2a(l-/3^)T 



(29) 



ifir) ^ (fa - {h- ahs) - J "m-zir) dr, 



where (/3(0) is the initial value of the angular variable and 



I3 = h+^^. 
a 



(30) 



Equations 129(1 and (|30l) represent a generalization of 
Eqs. (22) and (23) from Ref. 20] to the case of arbi- 
trary 6. As follows from Eq. H29|) . in the stationary case 
(r — > oo) and for vanishing easy-plane anisotropy, the 
system turns to the state with to^ = cos = ^ 1 in a 
complicated manner. At the initial time, the interac- 
tion energy has a jump dUa{hp — 0)/9t|^^p — ~2a(3h] 
further on, depending on the sign of (3, it increases or 
decreases in an oscillatory manner to reach the value of 
dU^{hp = 0)/5r|_^ = 0. 

Some analytical solutions can also be obtained for the 
stationary situations, when Eqs. (|BJ can be transformed 
into a trigonometric form with the general periodic solu- 
tions 6q and (po given by 



sin2iy9o = ^2/is//ip, 
cos^o = —h/Z{ipo). 



(31) 



The invariance conditions (O allow to select only two 
independent stationary states from those described by 
Eq. inH): 



^01 + = 7r/2, 



cos<^02 = -sign/is ^, 2 V" V " -■si-p 



l-^l-Ahl/hl), 



(32) 



sm (^02 

sin (poi = cos (^02 , cos (^oi = sm 'Pq2 , 
cos 6*01,2 = -ft./2'(v3oi,2)- 

As follows from Eqs. 132() . the stationary states are pos- 
sible only for hs satisfying the condition \hs\ < 0.5/ip. At 
the same time, the positive or negative value of the lon- 
gitudinal magnetization vector component is determined 
by the orientation of external magnetic field, being inde- 
pendent of the angle between the spin current and ez ■ It 
is important to note that the stationary states described 
by Eq. (|32|l are characterized by the potential energy of 
Eq. Isj, which is an integral of motion with the corre- 
sponding eigenvalues 



h'^/Z{ipQi 



(33) 



The differences between the equilibrium energies Uqs = 
1 + h and U04 — l + hp + h? /{i + hp) and these stationary 
states Usm,2 are 



Af/o 



llo<|h|<l+/i. 



Uoi — Uso 



il + hp){l + y + hpVl^)' 



AJ7o2|o<|/i|<i - Uso2 - Uo3 - y - Y 



(34) 



2hs 



V 



= ]-hp(l- Vl-a;2 



Therefore, Eqs. (|32|l to (|34l) demonstrate that neither sta- 
tionary states nor the corresponding energies depend on 
the damping coefficient a, as they are only determined by 
/i, hp and hs- From Eq. ((34|l follows that for hg = 0, the 
stationary state C/501 coincides with C/04, and Uso2 with 
U03. Under the influence of spin current, the state Uoi 
shifts towards the lower energy and becomes more ther- 
modynamically stable transforming into the state Usoi, 
while the state C/qs increases in energy and transforms 
into the state Uso2- As follows from Eq. (31), the rela- 
tive position of both stationary states on the energy scale 
could be then changed in controllable way by means of a 
proper choice of h and h^. 



D. Stability of the equilibrium and stationary 
states 

To investigate the stability of equilibrium (|10|l or sta- 
tionary (|32|l solutions, we subject them to small pertur- 
bations 

9^9, + 59,, ip = ipi + S(p^, (i = eorO) (35) 

with 69i <ti 9i, 6ipi <ti fi- Physically, such perturbations 
can have quite different origin: they can be caused either 
by a non-homogeneity of the sample, presence of con- 
tacts, fluctuations of the temperature or external fields;^ 
etc. Substituting (|35|l into © and performing a standard 
linearization procedure over small perturbations,"^'* one 
obtains the dynamic matrix a with the corresponding 
eigenvalues 

Ai,2 = 0.5 ^aii + a22 ± y/ (an — 022)^ + 4ai2a2i^ (36) 
and the matrix elements Uij 

ail — ^ cos 9i{a[Z(ipi) + h] + 0.5hp sin 2ipi + hs} 

+aZ{ipi) sin^ 
ai2 — hp sin. 9i{a COS 9iSin2(pi — cos2ipi), (37) 
021 = sin9iZ{ipi), 
0-22 = hp{cos 9i sin 2Lpi -I- a cos 2Lpi). 

After substituting the solutions H10|) into Eq. H37|) and 
then into Eq. (|36l) . one obtains three different solutions: 
a) for the state with sin^e = 0, (^e arbitrary: 

Ai = ^{a[Z{ipe) + h]+ 0.5hp sin 2ipe + hs}, (38) 
A2 = hp(zL sin 2ipe + a cos 2(pe), 
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with the upper and lower signs corresponding to 0^ = 
and 9e — tt, respectively. 

b) for the state with cos^e = —h, <Pe = 7r/2 



Ai,2 = ^ [/i + a(l + h - hp) 
±^[h hs + a{l + h + hp)Y + 4/ip(l - h^) 



(39) 



c) for the state with cos^?e — + hp), ipe = 0: 

0.5 



Ai, 



1 + k 



hhs + a{l + hp){l + 2hp + /i) 
±{[M, +a(l + /ip)(l + /t)]2 (40) 

-4/lp(l + /^p)[(l + /lp)2-/l2]}l/2l ^ 



The solutions of the characteristic equation, corre- 
sponding to the stationary states (|32|) . are rather com- 
plicated and we have calculated them numerically. It is 
worth noting that the solutions H38|) and (|39|l are real 
(Ai > and A2 < 0) - also for the stationary state (6'o2, 
9302)- In turn, the solutions (|40|l are complex-conjugated 
(Ai_2 = Ar ± iXirn), including the case of stationary state 
(6*01, foi)- 

The stability of the equilibrium or stationary states 
can be estimated using different stability criteria. In our 
case, we apply the method of Lyapunov^ evaluating the 
stability from the signs of Ai and A2, as well as the cri- 
teria of Gurvitz^^i^ taking into account the relations 
between the matrix elements Oy. In particular, accord- 
ing to Gurvitz, the sufficient conditions for stability of 
the system under consideration are 



an + 0,22 < 0, 



aiia22 > ai2a2i- 



(41) 



If Ai and A2 (see Eq. (13611 ) are known, the solutions 9{t) 
and (^(r) have the following form.'*'^ 
(a) for the case of real Ai and A2: 



- C'2ie^'^ 



Ci2e^-\ 
'C22e^^^.. 



(42) 



B{t) - B, 
v'yj) = 

(b) for the complex-conjugate Ai and A2: 

B{t) =d, + 2 [Dii cos(A„„t) - D12 sin(A„„T)] e^"", 
(^(r) = ipi + 2 [D21 cos(Ai„iT) - D22 sin(Aimr)] e^"-^ . (43) 

Here Cij and Dij are the integration constants, which are 
generally determined from the relevant boundary condi- 
tions. The latter ones can be written as 

^MIr=o = ^» + 27r, v{t)\^^^ = ip^ + ^r, 



for complex-conjugate Ai and A2. Equations 1421) and 
(I43|l can be rewritten as 



9i2 + 27rexp(A2T), 
= + 7rexp(A2r) 



(46) 



^2(r) = 
^2{r) 

for real Ai and A2, and 



6*1 (t) = e,i + 27rexp(A,,T) sin(A„„r), (47) 
= Vii +7rexp(A,.r) sin{\i„^T) 

for complex conjugated Ai and A2. As follows from 
Eqs. H46|l and (|47|l . the functions 0{t) and tf{T) can 
have either monotonic or oscillatory character, with the 
characteristics depending on the control parameters in a 
rather complex way, which can be investigated by numer- 
ical calculations described later. 



E. Phase averaging 

Let us come back to Eq. p3|l describing the time vari- 
ation of the interaction energy. Using Eqs. ©, one can 
rewrite this equation in the form 

2h,si,^'6^-^ = 2ash?e 

OT OT 

X {[Z{lp) COS e + hf + {0.5hpsm2ip + hs)^}. (48) 

On the right side of Eq. (|48() there is a positive quadratic 
form with a constant magnitude under any transforma- 
tion of the coordinate system. Hence, when the spin cur- 
rent is absent, the energy Uo{9, ip, r) decreases with time, 
i.e., the system tends to an equilibrium state. If h^ ^ 0, 
the variation of energy dUo / dr can be either positive or 
negative, depending on the sign and magnitude of the 
first term on the left side of Eq. (|48|l . As dtp/dr is deter- 
mined by hs for the fixed values of hp and a, it is natural 
that dU^/dT also depends mainly on the direction and 
magnitude of the spin current. 

To reveal the exact form of this dependence, one has 
to average Eq. (|48|l over the phase. When performing 
such an averaging, we take into account the periodicity 
of Eq. H48|l and also the invariance relations iQ. For the 
averaging procedure we will use two possible approaches: 
without any weighting coefficients^ 



(/) ^ (/(0,^)) - - / / f{0,v)ddd^, (49) 



de{T) 



dr 



= 2 



d(p{T) 



T = 

for real Ai and A2, and 



27r|A2| (44) g^j^j .^j^j^ weighting coefficients^*^: 



It=0 



d9{T) 



dr 



= 2 



difir) 



dr 



~ — 27r|A,;m,l 



Substituting into and ifSUI) . one obtains 
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a) for the averaging procedure (|49|l 



1 /dU^ 
a\ dr 



5^ 
32 ' 



(51) 



b) for the averaging procedure (|5 



4/i • hs 
3a 



dUo 
dr 



(1 



hi) 



(52) 

First of all, it is necessary to emphasize that the ex- 
pressions (|51|) and (|52|) have the same functional depen- 
dence on h and hs, correlating well with the invariance 
upon simultaneous replacement h —h and hs — > —hs 
(see Eq. 0). As one can see, the right hand sides of the 
obtained formulas do not include hs- This means that 
when the external magnetic field and the spin current 
have the same direction, the average energy can decrease 
faster than in the case of hs = 0. For antiparallel h and 
hs, the average {dUo/dr) can be negative, zero, or can 
become positive - but its absolute value is always smaller 
than \h ■ hs\ to obey Eqs. (|51|l or 1)52(1 . It is worth not- 
ing that the approximation 9 <^ I, used in (23 j leads to 
the violation of the symmetry condition (jSJ , which is also 
reflected in the averaging formula (|49|) for {dUo/dT): 



"^(l + ^ + h 
a \ 2 



1 



/dUo 

2oLdl \ dr 



(53) 



= 1 



hr. 



+ h] + 



The averaging procedure 1)50(1 in the given approximation 
should be changed to 



1 

7r?o 



9f{e,ip)d0dcp, 



(54) 







leading to two times greater coefficient at {dUo/dr) in 
(ESI- 

The difference between the expression ((53(1 and the cor- 
responding formula 1(16(1 from Ref. 20] is caused by the 
fact that we take into account the quadratic term in the 
expansion cosO\g^^ « 1 — 9^/2, which was neglected in 
Ref. 20]. 



III. RESULTS OF NUMERICAL 
CALCULATIONS AND THEIR DISCUSSION 

Using the formulas presented above, the magnetization 
dynamics of a single-domain ferromagnet can be investi- 
gated as a function of the control parameters: a, hp, hs 
and h. In our calculations, the number of control param- 
eters was reduced to h and hs by assuming a = 0.005 
(see Ref. [HI) and /ip = 5 (as in Ref. ^). 




-(1+hp) -1 1 +(1+hp) 

Magnetic field h [d.u.] 

FIG. 3: Field dependence of the stationary magnetic moment 
components. The curves labelled with 1, 2 and 3 correspond 
respectively to rrizi, rrixi and 40myi; whereas the curves 4, 5 
and 6 represent mz2, ^Omx2, and my2. (The curves 3 and 5 
have been multiplied by the factor of 40.) 



A. Analysis of equilibrium and stationary states 

Let us discuss first the dependence of stationary mag- 
netization states on the parameters h and hs- Fig- 
ure 3 shows the magnetic field dependence of nix — 
sin^ocosv^o, nT-y = sin sin c/?o , and = cos^o compo- 
nents of a unit magnetization vector m, corresponding to 
the stationary states l(32() . All other control parameters 
are kept constant in the calculations, namely a = 0.005, 
hp = 5 and hs = —0.03. As one can see, there is a range 
of magnetic field, where the magnetization components 
are non-zero for both stationary states. It is important 
to emphasize that the limiting field values dzhum depend 
on hg/hp for the stationary state (^oii V'oi), remaining 
practically independent on hs/hp for the state {9q2, ipo2) 
with hiirn = ±1. In both cases, the longitudinal compo- 
nent rUz varies linearly with magnetic field, whereas the 
transverse components nix < rUy have a non-linear field 
dependence. For better visual presentation of the data, 
the curves 3 and 5 in Fig. 3 have been multiplied by the 
factor of 40. 

The solutions of the characteristic equations for the 
equilibrium [(jSH), HOI)] and stationary [^^, (EHl, and 
(jSTjl ] states are plotted in Fig. 4 as a function of h. The 
solutions corresponding to the stationary state (6*01, foi) 
are complex conjugated (curves 1, 2 and 3), whereas the 
ones corresponding to (^02, f02) are real numbers of op- 
posite signs (curves 5 and 6). It is worth noting that the 
real part (curves 1 and 3) for the stationary state {9qi, 
foi) can be positive (area A), negative (area D) or equal 
to zero (point C). 

As concerns solutions of the characteristic equation for 
the equilibrium state ((39(1 . we note that the curves 5' 
and 6' in Fig. 4 practically coincide with the solutions 
for the stationary states {9o2, 'PQ2) (curves 5 and 6). For 
the state ((40(1 . only the imaginary parts of the solutions 
(curves 2' and 4') coincide with the imaginary parts de- 
scribing the state (6'oi , (poi ) (curves 2 and 4) , whereas the 
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Magnetic field h [d.u.] 

FIG. 4: Magnetic field dependence of tlie cliaracteristic equa- 
tion solutions: real (curves 1 and 3) and imaginary (curves 2 
and 4) parts of Ai, real parts of A2 (curves 5 and 6). Station- 
ary states: curves 1-6, equilibrium states: curves 1' - 6'. 
The curves 1,3 and 1', 3' are multiplied by the factor of 100. 
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FIG. 5: The real part of the characteristic equation solution 
as a function of h and hs'. 1 - stationary state, 2 - equilibrium 
state. 



real parts of these solutions differ significantly (curves 1', 
3' and 1, 3 in Fig. 4 are scaled by the factor of 100 for 
better visual presentation). It is important that for the 
given values of the control parameters, the real part of 
(|40|l does not depend on h. However, if the control pa- 
rameters are changed in a sufficiently wide range, the real 
parts of the solutions l|in)l and (Sqi, lpqi) behave similarly 
(Fig. 5). The difference between them is seen in the in- 
tersection of the surfaces \ri = /(/i, hs), {i = e, 0) with 
the Xr = plane, which for the stationary state takes 
place for weaker magnetic fields. 

According to the classification of specific points de- 
scribing the motion in the vicinity of equilibrium states?^ 
the solution (^02, ^02) corresponds to a saddle point 



FIG. 6: Bifurcation diagram of neutral stability for stationary 
(solid curves) and equilibrium (dashed curves) states for a = 
0.005 and different hp-. 1 and 1' correspond to ftp = 1; 2 and 
2' correspond to hp — 5; 3 and 3' correspond to hp = 10. 



(Fig. 4, B) with two attracting equilibrium states. The 
second stationary state (6101, ^poi) has a broader variety 
of possible motion types: unstable focus (Fig. 4, A), sta- 
ble focus (Fig. 4, D), and stable center (Fig. 4, C). It is 
interesting to note that the point, where the stable center 
appears is located on the boundary between the intervals 
of stable and unstable focuses. Therefore, the condition 
of reaching the stable center, Xr{hb) = 0, determines the 
bifurcation point hb (the boundary of neutral stability^) 
as 



T-bO 



where m = h'i — AKI . 

As follows from Eq. H55I) , with decreasing from zero 
to ha — ~hp/2, the value of /if, diminishes monotonously 

" =-(/ip/2a) + 
At 




from 



^{l + hp){l + 2hp) to /i^™ 



^(/ip/2a)2-t-(l-|-0.5/ip)2 (see Fig. 6). At the same 
time, /if,o shows a non-linear increase with the anisotropy 
hp. It is worth noting that for the equilibrium state H4U|I . 
the neutral stability field boundary h^e (Fig. 6, dashed 
lines) is determined from the expression. 



(l + fep)(l + 2fep) 

h^a^^ + I + h„ 



(56) 



satisfying also the condition h^e > hbo- In this way, l|55|l 
and (|56|l determine the system stability limits, as the 
corresponding states are stable when h > hbe (i — 0, e) 
and unstable in the opposite case. 

The two-dimensional plot of the energy of both Usoi 
and Us02 stationary states, calculated using H34|l . is pre- 
sented in Fig. 7 for < \h\ < 1 (state (6*02, ¥^02)) and 
for <\h\ <6 (state {601, Voi)); < a; = \2hs/hp\ < 1. 
For h = both surfaces merge, taking the energy value 
Uso{h = 0) = 1 -I- 0.5/ip. With increasing h, the energy 
of the low-energy state Us02 increases slightly, but this 
increase is faster for larger values of the spin current h^. 
Contrary, the changes of the high-energy state Usoi with 
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FIG. 7: Dependence of the stationary energy on h and x 
\2hs/hp\: 1 - ;7soi, 2 - Us02. 
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FIG. 9: Dependence of mz on external magnetic field for the 
equilibrium state at different observation times. 
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FIG. 8: Dependence of on external magnetic field for the 
stationary state at different observation times. 



increase of either h and hg are significant, non- linear, and 
lead to much greater difference in the energy of the states 
Usoi and Us02- Therefore, we will focus below on the non- 
stationary magnetization dynamics, corresponding to the 
excited state Usqi- 

The time and field behavior of the longitudinal mag- 
netization component rriz under small perturbations are 
illustrated in Figs. 8 and 9 for the stationary (^oii foi) 
and equilibrium {U04) states, using formulas (|^ . l|H7|) 
and H47|) . At the initial moment (r — 0) and after a 
short period of time (r = 1), the changes of the com- 
ponent with the increasing magnetic field are equal for 
both stationary and equilibrium states, but with a fur- 
ther increase of the observation time, the character of 
these curves becomes different, at first slightly (t = 10) 
and then significantly (t > 10^). ft is worth noting that 
the oscillations become chaotic, causing niz to take ran- 
dom values in the range of — 1 < tUz < 1. It is important 
that for T > 10^ and h > 4.2 (Fig. 8), the perturbation 
of the state (6*01, V^oi) relaxes to the linear = f{h) de- 



pendence, comparable with Fig. 3 (curve 1), while the re- 
laxation of equilibrium state 1)401) takes place for r > 10^. 

Therefore, under the influence of small perturbations, 
the high-energy stationary and equilibrium states of the 
ferromagnetic system in magnetic field are characterized 
with the complex temporal dynamics, which can be in- 
vestigated in detail using the numerical methods to solve 
Eqs. ®. 



B. Results of numerical modelling 

The type of spin orientation in the system varies signif- 
icantly with a change of the applied magnetic field. Our 
investigation reveals the presence of several oscillation 
regimes depending on the value of the applied magnetic 
field h. The set of differential equations © was solved 
numerically using the Runge-Kutta method of the fourth 
order. This allows to observe the motion of the magne- 
tization vector in the three-dimensional space, forming a 
phase portrait of the system. To determine the effect of 
magnetic field on switching between the states rriz = ±1, 
the magnetization vector was assumed to point down at 
the initial moment, i.e., = my = 0, = —1. 

The most characteristic phase portraits of the sys- 
tem, shown in Fig. 10, were obtained for hg — —0.03, 
— 5 < ft. < 15, a = 0.005, and hp = 5. The density plot of 
ruz component versus applied field (Fig. 10(a)) features 
darker areas corresponding to the most frequent val- 
ues at a given magnetic field. The field dependence of 
the power spectral density S{u!) of the component 
is given in Fig. 10(b), showing the peaks as dark areas. 
The precession velocity can be estimated from the tra- 
jectory tracing curve evolutionj^SiSi shown in Fig. 10(c). 
Here the dark areas correspond to the moments of time r, 
when the vector connecting two consecutive phase points 
becomes parallel to the one between two initial points; 
the distance between two black stripes At is inversely 
proportional to the velocity of the phase point. The 
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FIG. 10: Phase portraits and main system characteris- 
tics vs apphed magnetic field: (a) rriz density plot, (b) 
power spectral density S{ui), (c) evolution of trajectory trac- 
ing curve, and (d) Hausdorff dimension Dh and maximum 
Lyapunov exponent Xmax- 



Hausdorff dimension Dh of the phase portrait and max- 
imum Lyapunov exponent \max are shown in Fig. 10(d) 
for different values of the apphed field. 

For magnetic fields h < 1 the phase point remains in 
the vicinity of the ground state = — 1 during all the 
observation time, which is clearly seen at the density plot 
as the dark horizontal bottom line. The trajectory trac- 
ing curves allow to discern a negligibly small movement of 
the phase point characterized by a period increasing with 
the applied magnetic field. The tendency to instability is 
also notable from the \max curve reaching positive val- 
ues. When the apphed field h surpasses h — I, the system 
switches to the magnetization precession mode, charac- 
terized by the periodic phase point orbit in mx^ vriz plane 
(Fig. 10, h — 1). The corresponding place of the density 
plot shows two darker branches, marking the upper and 



lower limits of values, which are contracting towards 
Tfiz = — 1 with increasing magnetic field. Darkening of 
almost homogeneous gray area between these branches at 
higher h is caused by the limit cycle stability loss, visu- 
ally observable as a phase trajectory broadening (Fig. 10, 
h—\\ At a certain merging point of /i w 5, the darker 
branches join together and the phase portrait of the sys- 
tem turns to the focus, i.e., the magnetization precession 
converges to a value close to the state TOz = — 1 . 

It is worth noting that the S(io) plot reveals several 
peaks corresponding to the harmonic oscillations in the 
limit cycle mode. The number of these peaks dimin- 
ishes upon reaching the magnetic fields above which the 
phase portrait of the system turns into a focus. As can 
be deduced from the trajectory tracing plot, the initial 
shrinking of the limit cycle oscillation mode is accompa- 
nied by the precession frequency increase almost up to 
the merging point. The movement of the phase point 
becomes then slower for higher magnetic fields. The 
limit cycle stability loss and the corresponding increase 
of the phase portrait density could be clearly seen from 
the Hausdorff dimension curve, featuring a rapid increase 
from 1.55 to 2.05 before the merging point and almost no 
changes above it. The Lyapunov exponent has a noisy 
maximum at /i « 3, revealing a return to the stable state 
upon switching to the converging magnetization preces- 
sion mode. 

It is worth noting that in the whole range of applied 
field, 6 < /i < 11, the magnetization vector is restricted 
to the closest vicinity of its ground state = — 1. At 
the same time, as it can be deduced from the trajectory 
tracing plot, the frequency of the phase point precession 
increases, and starting from h k, 11.1 the system begins 
to tend to the upper state with to^ — 1 through the mag- 
netization precession. The phase portrait corresponding 
to this oscillation mode represents a sphere formed by a 
tight spiral with Hausdorff density Dh — 2.3 (Fig. 10, 
/i=14). It is worth noting a rather stable nature of such 
precession, described with almost zero values of Xmax ■ 

A similar magnetization dynamics is also observed for 
different values of /ig, as it can be seen from Fig. 11, pre- 
senting several density plots. In all the cases, the 
system under consideration does not leave the vicinity of 
the ground state to^ = — 1 until the magnetic field over- 
comes the value h = \. Above this point, it is possible 
to observe a magnetization vector precession to the state 
mz = 1 for hs = and hs — —0.01, while for higher neg- 
ative spin current values, the system is switched to the 
limit cycle mode. From now on, the system evolution 
follows the same scenario with the limit cycle contract- 
ing into a focus. Depending on the spin current value, 
the system remains in the vicinity of the bottom ground 
state for different ranges of the applied field h. Starting 
from some threshold field, which increases with the ab- 
solute value of hs, one can observe the switching of the 
system from the bottom to the upper ground state. 

In the framework of our description, we have performed 
the averaging of the expression H48|) over time according 
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FIG. 11: The system evolution with the apphed magnetic 
field h: rriz density plots for different values of the spin cur- 
rent. 
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FIG. 12: Dependence of a~^{dUo/dT) on the applied mag- 
netic field for different averaging procedures: 1 - according to 
JETl . 2 - from the formula lICT . 3 - using lIH^ . 4 - from JETl . 
5 - according to the formula 15311 . 



to Ref. m 



/dUo\ 




= lim — ; 


\dr / 


T^^T J 



dr 



dr, 



(57) 



where the value dUo/dr was calculated for 6{t) and (/^(r), 
determined from the numerical solution of Eqs. © using 
the Runge-Kutta method, and changing the condition 
T — > cx) in (|57|) to T — > T^axy where T^ax was selected 
from {dUo /dr) saturation condition. 
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FIG. 13: Time evolution of (curves 1, 2, 3) and dUt^/dr 
(curves 1', 2', 3'), calculated for hs = —0.03 and different 
applied magnetic fields: 1,1' for h = 11; 2,2' for h — 13; 3,3' 
for h = 15. 
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FIG. 14: Magnetization component (curves 1-4) and 
dUo/dr (curves l'-4') as a function of magnetic field for 
r — 800 and different spin currents: 1,1' - for hs = —0.05; 2,2' 
- for hs = -0.04; 3,3' - for hs = -0.03; 4,4' - for hs = -0.02. 



^From the comparison of the curve 1 in Fig. 12, calcu- 
lated according to formula (|57|) . with the curves 2 and 3 
calculated from (|51|l and (|52() . follows that there are sig- 
nificant quantitative and qualitative differences between 
the time (|57|l and phase ll^ . (|5n|l averaging procedures. 
However, for the Sun's model (9^1) with h > 0, the 
results of phase averaging (Fig. 12, curve 5) calculated ac- 
cording to the formula H53|) approach asymptotically the 
results of time averaging (curve 4, obtained for 6*0 = 1) 
with increasing external magnetic field. Such a differ- 
ence between the phase and time averaging of Eq. ifT^ 
for the case of arbitrary 9 and for the particular case of 
9 <^ 1 are caused by the fact that the averaging pro- 
cedures (|49f) and (|50|l consider incorrectly the weighting 
factors for the polar angle 9. One of the other evidence 
supporting this point is the expression for charge current 
Ic, determined on the base of (|51|) (or H52|l 'l from the 
condition {dUo/dr) = 0^ 
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FIG. 15: Dependence of switching time At on /i for different 
spin currents: 1 - for hs — —0.030, 2 - for hs = —0.035, 3 - for 
hs = —0.040. Behavior of the parameters Atoo and hr with 
spin current hs is given in the inset. 
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which reveals an unexplainable divergence at /i — > 0. At 
the same time, expression ()54(l practically corresponds to 
the averaging only over the angle ip, which - taking into 
account a good correlation of curves 4 and 5 (Fig. 12) - 
seems to be the most appropriate one. This result allows 
to simplify significantly initial equations lOl , © j and Ill3|l 
by averaging their right-hand sides over ip. As a result, 
instead of ^ one obtains 



dr 



dr 



= a{l + 0.5ft,p)(l - ml){m 

dip 
97 ^ 



(1 + 0.5/ip)TO^ 



ah, — h. 



(59) 



Instead of Isj and ((T^ it is possible to get correspond- 
ingly 



C/o = (1 + 0.5/ip)(l - ml) - 2hm,, 



(60) 



and 

dUo 
dr 



-2a{l + Q.5hpf{l-ml){m,+/3p){m,+^), (61) 



with niz = cos9, f3p = {h + ha/a)/{l + 0.5/ij,) and 7 = 
/i/(l + 0.5/ip). 

The ordinary differential equations (|59|) are of 
separable-variable type and allow to find solutions for 
the averaged azimuthal angle Tp at known mz{T) in the 
form 



ip{T) = Tp{Q)-{h-ahs)T- ( 1 + ^/if 



m,(T)dT, (62) 



where ^(0) is the initial value of the angle (i.e., its phase). 

In the particular case of (3p — ±1, mzir) can be found 
from the following transcendent equation: 



arctanh niz T 



1 



1 ± TUz 



(63) 



arctanh rn^o T 



1 



1 ± m^o 



T2a (1 + 0.5/ip) r. 



where nizo is the initial value of the longitudinal magne- 
tization vector component. As follows from Eq. H63|l . the 
sign change of the parameter /3p is equivalent to the si- 
multaneous sign change of and r. It is worth noting, 
that despite the case /3p = ±1 was analyzed in detail2£ 
(see, for example, Eq. (19)) using numerical methods, the 
analytic equation for mzij) was not presented there. 

For (ip 7^ ±1, rnzir) satisfies the following transcen- 
dent equation: 



0P 



{1 ~ ml){mz + I3p)- 



1 — n 

1 + rUzQ 
1 - mzo, 
X exp [-2a(l -I- 0.5/ip)(l - /3^)t] 



(l-"^zo)("^^o + /5p)" 



(64) 



It is easy to notice that expressions (|62|l and H64|l are 
generalizations of formula (|29() for the case of hp ^ 0. 

As the case /3p = ±1 (Eq. jHSJ) was already investi- 
gated in detailfSS where it was used for determination of 



the threshold field 



,(±) 



for the switching from parallel 

to antiparallel ruz state, we will consider below only the 
switching peculiarities defined by Eq. H64|l . 

In Figs. 13 and 14 we present the time and field de- 
pendence of niz and dUo/dr, calculated according to 
Eqs. H64() and (|61|l . respectively. The dashed areas in 
Fig. 13 correspond to mz{T) values obtained by the nu- 
merical integration of the system using the Runge-Kutta 
method. As one can see, for different values of h, mz{T) 
determined from Eq. H64|l forms the outline of numer- 
ically calculated magnetization vector projection. It is 
worth noting that the good correlation between the an- 
alytical and numerical results can be also observed for 
the function mz{h) for different values of spin current hs 
(Fig. 14). The latter proves that the proposed simplifi- 
cation of initial equations 10 , © , and ifT^ by averaging 
over ip is correct. The dashed curves in Figs. 13 and 14 
correspond to time- and magnetic field-dependent val- 
ues of dUo/dr, calculated on the base of Eqs. and 
(|64|l . Our analysis has shown that these curves can be 
approximated well with the Gaussian distribution.^-'^ The 
values mzmih, hs, hp, a) corresponding to the minimum 
of dUo/dr, can be found by solving the cubic equation 



1 



+ -(/3p + 7)m^„ + -(7/?p-l) 



.-^(/3p+7) = 0, 
(65) 

which can have one or three real solutions, of which only 
those satisfying the condition — 1 < mzm < 1 should be 
taken into account. 

As follows from Fig. 13, at the fixed values of hs, hp 
and a, the minimum value of dUo/dr = f{h) decreases 
exponentially with time, while the similar minimum of 
dUo/dr = Tpihs) at fixed t, hp, and a grows linearly with 



14 



the increase of hg (Fig. 14). It is important to empha- 
size that in the latter case the value of applied magnetic 
field, at which turns to zero, decreases linearly with 
increasing hs- 

As one can see from Figs. 13 and 14, the switching 
time between the ground states ruz = — 1 and niz = +1 
has the same order of magnitude as the half-width of 
the dUo/dr peak. This half-width remains practically 
unchanged with h at the fixed r, hp, and a (Fig. 14), but 
depends significantly on the magnetic field with fixed hg , 
hp, and a (Fig. 13). As it turns out (Fig. 15), the half- 
width AT{h) of the peak dUo/dr = f{h,T) diminishes 
exponentially with increasing h according to formula 

At = Atoo + (Aro - Ar^o) exp{~h/hr), (66) 

where the parameters Atoo, Atq and hr depend on hs, 
hp, and a. It is worth noting that for the fixed applied 
magnetic field h, the increase of kg leads to decrease of 
the switching time At. The dependence of the param- 
eters Atoo and h,. on the spin current kg are presented 
in the inset to Fig. 15. As one can see from this figure, 
Atoo and hr are characterized by a linear decrease and 
increase with growing hg, respectively. The results pre- 
sented in Fig. 15 show that the increase of the absolute 
value of spin current leads to destabilization of the sys- 
tem, which now needs more time to switch between the 
ground states niz — —1 and = +1. To the contrary, 
the applied magnetic field has the stabilizing action, lead- 
ing to a decrease of Ar. It is also worth to emphasize that 
for h > 25, the infiuence of the applied magnetic field be- 
comes dominant compared to the role of the spin current 
hg. It is worth noting that formulas (|64|l and H65|) allow 
to obtain analytical expressions for the parameters (|66|l . 
making it possible to investigate in detail the dependence 
of the switching time Ar on h, hg,hp, and a. 

IV. CONCLUSION 

We have studied the magnetization dynamics of a fer- 
romagnetic system subject to the spin-polarized current. 
We have used the methods of non-equilibrium thermody- 



namics, which have been developed to describe the self- 
organization processesi^~i^£'Sfi'^ Our results show that 
the ferromagnetic system displays a complex dynamics 
of instabilities with applied magnetic field and spin cur- 
rent. We have demonstrated that the method can be 
used to describe the dynamical properties of ferromag- 
netic nanostructures important for spintronics applica- 
tions. 

The behavior of equilibrium and stationary states has 
been investigated for a wide range of external parame- 
ters, without any limits on the angular variables describ- 
ing the system magnetization. It has been shown that 
under certain conditions, the system can be switched to 
the oscillation mode regime with negligibly small damp- 
ing. The phase portrait evolution has been investigated 
in detail for different values of external magnetic fields. 
The obtained results demonstrate a possibility to control 
the operating modes of the ferromagnetic components in 
spintronic devices by a proper choice of the ferromagnetic 
material with appropriate easy-plane anisotropy and the 
damping coefficient. It has been shown that the averag- 
ing of the initial system over the azimuthal angle allows 
to obtain analytical expressions, which make possible the 
detailed investigation of the switching peculiarities of the 
longitudinal magnetization between the states = — 1 
and rriz = -|-1. The switching time has been analyzed 
in dependence on the external magnetic field and spin 
current. 

The obtained results proves the successful applica- 
tion of self-organization methodology to solve theoretical 
problems of spintronics, which can yield new quantitative 
and qualitative results. 
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